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Abstract 

The Aircraft Morphing Program at NASA Langley 
Research Center explores opportunities to improve 
airframe designs with smart technologies. Two 
elements of this basic research program are 
multidisciplinary design optimization (MDO) and 
advanced flow control. This paper describes examples 
where MDO techniques such as sensitivity analysis, 
automatic differentiation, and genetic algorithms 
contribute to the design of novel control systems. In 
the test case, the design and use of distributed shape- 
change devices to provide low-rate maneuvering 
capability tor a tailless aircraft is considered. The 
ability of MDO to add value to control system 
development is illustrated using results from several 
years of research funded by the Aircraft Morphing 
Program. 

Introduction 

Researchers who specialize in multidisciplinary design 
optimization (MDO) and those who specialize in 
optimal control have a natural affinity, since all use 
mathematical optimization techniques. Yet. few 
examples ol MDO research including flight control as 
one of the disciplines exist. This apparent 
contradiction stems from the traditional process of 
designing aircraft night control systems. Traditional 


aircraft conceptual design (see for example Ref. 1) 
assumes a small number of conventional control 
devices: ailerons (roll effectors), elevator (pitch 
effector), and rudder (yaw effector). The size and 
location ot these devices can be estimated by using 
historical databases for component weight and control 
effectiveness information. The detailed control law 
design is postponed until the aircraft configuration is 
frozen and until precise control moments are measured 
or predicted. 

The Aircraft Morphing Program envisions new types of 
control devices, such as inflatable bladders or 
oscillatory jets distributed over the wing surface, which 
achieve control by changing the real or virtual shape of 
the wing. - These novel control systems require an 
equally revolutionary control design process. This 
paper suggests that MDO techniques, such as 
automatic differentiation for calculating sensitivities 
and genetic algorithm (GA) optimization procedures, 
are essential components of that new process. 

Control System Design Process 

The proposed control system design process is 
illustrated by application to a tailless fighter aircraft 
concept. The numerous control devices (effector 
arrays) were modeled as generic shape-change devices 
that respond immediately to commands from the 
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controller. This kind of simplification was justified 
because the aim of the paper is to demonstrate MIX) 
techniques. 

Figure l illustrates the proposed control system design 
process. In Fig. 1 , the arrows indicate the flow of data 
and the boxes represent steps in the design process. 
For example, the first step is to develop the vehicle 
concept and the final step is to develop and test the 
control law. The arrow which connects the Final box to 
the initial box indicates that unsatisfactory results from 
the control law design may necessitate changes in the 
vehicle concept. In Fig. 1, the boxes with a dark border 
represent steps where MDO techniques can add value. 
The boxes with a light border represent steps where 
engineering judgement is especially important. Each 
step in the process is described in this section. 

Development of Ve hicle Concept 

The first step in this multidisciplinary control system 
design process was to create a computational fluid 
dynamics (CFD) model of the vehicle and the control 
devices. In the present research, generic shape-change 
devices controlled a representative aircraft 
configuration called ICE (Innovative Control 
Effectors), created by Lockheed Martin*. The ICE 
design, described in Refs. 3 and 4 and shown in Fig. 2, 
was used under a cooperative agreement with 
Lockheed Martin. The effectors were shape-change 
dev ices modeled as bumps on the surface of the w ing. 
For the current proof-of-concept studies, existing grid 
points were deformed in the direction of the surface 
normal to represent potential shape-change devices. 
The CFD model of the ICE configuration is evaluated 
by an aerodynamic panel code called PMARC. 

Prediction of Control Moments 

The second step was to predict the sensitivities of the 
control moments to a shape change at each grid point 
of the CFD model. The sensitivities were obtained by 
automatic differentiation ot PMARC with the ADIFOR 
code generation tool/ 1 The sensitivities required were 

t h e partial derivatives of roll (C,), 

\ dh dh dh ) 

pitch (C,„), and yaw (C„) moments with respect to a 
displacement /?. The details ot calculating a surface 
normal at every grid point and applying a change in 
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height h along that normal are described by Park in 
Ref. 7. For the present ICE configuration model, these 
calculations required the derivatives of 3 output 
quantities with respect to 1394 input quantities. 
Calculating such a large number ol derivatives was 
practical because of the adjoint option in version 3.0 of 
the ADIFOR software. 8 Park estimates that he spent a 
total of about one week modifying the PMARC code 
and applying ADIFOR-3.0. Calculating the 
sensitivities required about one hour ol CPU time on a 
high speed engineering workstation. 

Given a 3x1394 matrix of partial derivatives for grid 
locations on the right wing, the partial derivatives for 
corresponding locations on the left wing were 
constructed by assuming they had same magnitudes but 
the roll and yaw derivatives had opposite sign. Given 
this matrix of all partial derivatives, the moments m 
resulting from activation of any set of effectors was 
estimated with a matrix multiplication: 

m=[B]A cmU (l) 

where u tm<l is an //-vector of device heights h and the 
matrix B is constructed by selecting the n columns 
associated w ith those dev ices. Note that this method of 
estimating moments implies linear superposition, 
which neglects control effector interactions. On the 
other hand, if m nHlt represents a vector of desired roll, 
pitch, and yaw moments, then w r/Mf /Can be calculated by 
using the pseudo-inverse allocation method discussed 
in Ref. 9: 

< 21 

where T denotes the matrix transpose. As in Ref. 9. 
heights were restricted to positive values less than 
some maximum achievable device height. 1 herefore, 
if any element of u (Wtl was negative that element could 
be set to zero, and the corresponding element on the 
opposite wing could be increased by the same value. 
Then the resulting vector was normalized so that the 
maximum element equaled the maximum achievable 
device height. Finally, Eq. (1) was evaluated and the 
target moments ni tW d were compared w ith the 
achievable moments m. 

Definition of Effector Array Candidates 

Equations (1) and (2) were used to estimate the control 
moments produced by an array of effectors on the K E 
vehicle. An interactive design tool was created to let a 
researcher quickly build up and analyze potential 
locations (i. e., effector arrays) by selecting grid points 
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and assigning device heights. This interactive tool, 
which uses MATLAB software developed by The 
MathWorks, Inc., helps the control law designer to 
determine good potential locations so that each device 
can produce the forces and moments required to 
maneuver the vehicle. Figure 3 illustrates the graphical 
user interface for the MATLAB -based effector array 
design tool. The wing planform on the left side of the 
display presents the CFD grid. The three contour plots 
on the right side of the display present the sensitivity 
data produced by step 2. The top contour plot indicates 
roll moments on the top and bottom wing surfaces, the 
middle plot indicates pitching moments, and the 
bottom plot indicates yaw moments. 

Once an array of effectors were defined by specifying 
locations and heights, the designer obtained a 
preliminary prediction of the arrays effectiveness based 
on Eqs. ( I ) and (2) and ADI FOR sensitivities. Next, he 
generated a perturbed geometry grid that included the 
deployed effector array. This geometry file could then 
be reevaluated with the PMARC aerodynamic analysis 
program. New sensitivity derivatives were calculated 
based on the perturbed PMARC’ results. These 
PMARC-based sensitiv ities could be used with Fqs. ( 1) 
and (2) to further assess the effectiveness of the array. 
For the present study, perturbed dev ice heights were 
less than or equal to 0.2 feet. 

Figure 4 displays a 17x41 grid of candidate device 
locations on the right upper wing. Circles indicate one 
selected effector array. This array contains 28 dev ices. 
A right and left pair of these wing tip arrays can 
provide the estimated roil, pitch, and yaw moments 
shown in the first column of Table 1. Reanalysis of 
these same devices with a deformed grid in PMARC’ 
gives slightly modified estimates shown in the second 
column of Table 1. In this way, 34 different effector 
arrays were selected and evaluated. The two methods 
of estimating moments were quite consistent; in 
general the ADIFOR estimates were smaller than the 
PMARC estimates. Therefore, reliance on ADIFOR 
estimates for optimization should produce a 
conservative design. 


Table I. Control Moment Estimates Compared with 
PMARC’ Reanalysis Results. 


Moment 

ADIFOR 

PMARC 

Roll 

-L44E-04 

-1.56E-04 

Pitch 

-2.02E-04 

-2.20E-04 

Yaw 

-0.22E-04 

-0.28E-04 


The PMARC estimates for the effector arrays are 
plotted in Fig. 5. Estimates based on Eqs. (1) and (2), 
with sensitivity information generated by ADIFOR, 


would create a very similar plot. The candidate 
effector arrays were located in 7 regions: on the upper 
and lower leading edge (LH), on upper and lower 
trailing edge (TE), on upper and lower wing tips, and 
on the upper surface near the middle of the wing. The 
arrays were not necessarily disjoint. In fact, some 
devices were members of sev eral effector arrays. 

The ellipses in Fig. 5 indicate a suite of four effector 
arrays that were studied in Ref. 9. This suite was used 
in a six-degree-of-freedom dynamic simulation to 
investigate the unaugmented and augmented aircraft 
dynamics. Results of that simulation are reported in 
Ref. 9 and indicate that a 1 0-degree- per-second roll 
rate is the maximum achievable with this suite of 
effectors. Those results suggest that this particular 
suite could be valuable for mild maneuvering or could 
be used in an autopilot to keep the wings level but 
could not take the place of conventional control 
surfaces. 

Selection of Optimal Effector Arrays 

The manually selected effector suite studied in Ref. 9 
required 82 individual devices and did not completely 
meet the goals set by the designer. Given enough time 
and good intuition, the designer might have selected 
other effector suites with better characteristics. 
Alternately, the initial exploration for candidate 
effector suites can be accomplished with discrete 
optimization techniques. For example, Ref. 10 
contains a literature survey of actuator placement 
research that indicates good results for a wide range of 
applications. 

With the MATLAB-based tool, the designer can define 
a large number of arrays and then can predict their 
control moments by using the PMARC analysis code. 
However, selecting the best set of arrays from this 
potential pool is a combinatorial problem that can be 
solved by using MDO techniques. The current study 
selected a GA-based optimization approach. The goal 
of the optimization was to reduce the number of 
devices required and to satisfy control effectiveness 
criteria. The GA is not suggested as a replacement for 
the designer, but as a tool for screening potential 
effector arrays and thus allowing the designer to 
consider only the most promising. 

A multilevel GA, used to select control device 
locations, is described in Ref. I 1. The goal of the GA 
is to select the minimum set of dev ices that can provide 
the required uncoupled control moments (e. g., provide 
sufficient pitching moment without adverse roll or 
yaw). That GA technique was tuned and validated 
with a simplified wing model." 12 
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In the current study, the GA described in Kef. 1 1 was 
adapted to select effector arrays on the ICE model. 
Again, the objective was to find the minimum number 
of devices required to provide uncoupled roll, pitch, 
and yaw moments. Each member in the GA population 
represented one possible effector suite. Each 
individual was evaluated three times to determine roll, 
pitch, and yaw moments and to compare these 
achievable values with the target values. Failure to 
meet any of these targets caused a penalty to be added 
to the fitness (i. e., objective) function. 

In Ref. 11, the penalty has a fixed size. After some 
experimentation, a step-linear penalty was determined 
to work better. Thus the fitness function, J, can be 
written as follows: 


this problem, the GA was allowed to pick at most one 
array from each of the 7 regions shown in Fig. 5. The 
GA implementation would be simplified if each region 
had the same number of arrays. So, because the upper 
TE region contained 8 arrays, arrays in the other 
regions were duplicated until each region contained 8 
arrays. In this way, the string length for the phase I 
GA was set to 7 and each digit in the string could have 
a value between 0 and 8. Thus, the number of possible 
combinations was 9 7 (approximately 4,800,000), 
although not all of these combinations represented 
unique effector suites. The duplication of effector 
arrays should not have had a significant impact on the 
convergence of the GA. Although, duplication does 
produce many members of the population with the 
same fitness value, a GA is especially well suited for 
optimization problems with this characteristic. 


J = n + 



+ w\ 


o 

c„ 


+ w 



( 3 ) 


where n is the number of devices, w, are the minimum 
penalties, C, are the moments and C,* are the targets. 
For cases reported in this paper, w, = 150 if the 
moments are less than the targets (e. g., if C, < C f * ) 
and w, = 0 otherwise. The number 150 is an 
appropriate size for the minimum penalty because it 
has the same order of magnitude as n judging from the 
number of devices in the manually-selected suite. 
Thus, the GA will be encouraged to drive all the 
penalty terms to zero by choosing a sufficient number 
of devices. 


For the current study, the GA was developed and tested 
in two phases. During phase I, the 34 effector arrays 
selected by Raney were defined as the set of possible 
arrays. The GA was allowed to select up to 7 arrays 
from the original 34. The effectiveness estimate for 
each suite of effector arrays was based on the PMARC 
reanalysis data used in Ref. 9. The target values (see 
Eq. 3) were set to C f * = 6.0E-04, C„* = 5.0E-04, and 
C„* = 3.0E-04. During phase II, a set of 349 
individual devices could be selected independently. 
The effectiveness estimate for each suite of devices 
was based on the ADI FOR data. The target values for 
phase II are set to C,* = 5.0E-04, C,„* = 5.0E-04, and 
C„* = I.0E-04; the targets were reduced because the 
estimates based on ADIFOR data consistently 
underestimated roll and yaw. 


The GA population size for phase I was 200. An initial 
population of members was produced randomly and 
their fitness evaluated. Successive generations were 
produced by the GA operations of tournament 
selection, uniform crossover, and mutation, with a 
mutation rate of 5%. The maximum number of 
generations was set to 300. A single execution 
consisting of 300 generations of the GA procedure 
requires about one hour on a engineering workstation. 
For complete descriptions of GA techniques and 
definitions of GA parameters, see Ref. 13. 

Typical results of the GA are shown in Fig. 6. Notice 
that the scales on each figure are different in order to 
emphasize several points about the convergence 
history. Figure 6a shows the value of the fitness 
function averaged over all population members in each 
generation. The maximum and minimum fitness is also 
plotted. Figure 6a shows clearly that the population 
maintained adequate diversity for 300 generations. 
Figure 6b shows the average and minimum fitness for 
the first 20 generations. Notice that the average fitness 
reduced dramatically over this interval and approached 
the value of 150. This trend suggests that all members 
of the population were converging towards designs 
whose achievable moments were close to or better than 
the target values. Figure 6c shows the best fitness ever 
calculated as a function of generation number. Figure 
6c indicates the original population contained at least 
one member with a fitness less than 150, suggesting 
that some members met all the performance targets. 
Notice also that the best individual had 96 devices and 
was found before generation 20. 


Obviously, phase I was a much smaller combinatorial 
problem, but it had its own complications. For 
example, some of the arrays overlapped and so not all 
possible combinations were allowable. To circumvent 


Even though 300 generations and about 60,000 
different effector suites (i. e., about 1.25% of all 
possible combinations) were evaluated by the GA, 
apparently no combination w ith fewer than 96 devices 
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could be found. This global best design is pictured in 
Fig. 7. The design had three arrays on the upper wing 
surface and two on the lower surface for a total of 96 
devices in 5 arrays as compared to 82 devices in 4 
arrays selected by the manual method. However, this 
96-device design met all the targets while those 
selected by the manual method did not. 

Based on encouraging results from phase I, the GA was 
used to select individual devices from a set of 349. 
This set contained all of the devices that made up the 
34 arrays in phase I, plus some other devices that 
seemed promising. 

The phase II GA was tested with a variety of crossover 
strategies and mutation rates. A typical execution of 
the phase II GA requires about 8 hours of CPU time. 
The choice of strategy does not seem to affect the 
results or the efficiency very much. The results in Fig. 
8 were produced with a single -point crossover and a 
mutation rate of I %. 


The members of the phase II GA population were 
initially assigned random numbers from 0 to 349. If 
the same random number was generated more than 
once in any member, then all the duplicates were set to 
zero. Duplicates were similarly set to zero following 
crossover and mutation operations. The string length 
was set to 100, so the number of possible combinations 


of 100 devices chosen from 349 was 



4 x I0 tw . 


Because the number of combinations was much larger 
than in phase I, the population size was increased to 
300 and the maximum number of generations was set 
to 500. Thus, about 150,000 individual members were 
evaluated during each repetition of the GA. This 
number represented a tiny percentage of all possible 
combinations. 


Figure 8 shows typical convergence performance for 
the GA, and Fig. 9 shows the 45 devices selected by 
this execution of the GA. Notice in Fig. 8, that even 
the ’’best ever” fitness value was initially above 500; 
this high initial value means that all individuals in the 
population were heavily penalized. After about 100 
generations, the best fitness dropped sharply, indicating 
that a design meeting the targets had been found. 


proportional fashion: each dev ice in an effector array 
was set to the same height. As greater moments were 
required by the control system, the height of a given 
effector array was increased until the limiting height 
was reached. By using these arrays, the control system 
was able to stabilize and maneuver the vehicle without 
conventional moving surfaces such as ailerons or a 
rudder. The predicted authority of these devices was 
still rather low when compared with that of a rudder or 
aileron, so the control system generated relatively low- 
rate maneuvers (roll rates of 5 to 10 degrees per 
second). Future research will focus on experimental 
validation of the predicted authority of various flow 
control devices and on better estimates of their 
effectiveness. 

Fi gure 10 contains time history plots comparing the 
phase I GA suite (thick solid), the original manually- 
selected suite (thin dashed), and the ideal moment 
targets (thin solid). The roll and yaw moments coming 
from the GA suite look very good-slightly better than 
those from the original suite. Both suites cause an 
undesired pitch perturbation (see Fig. 1 0b), but the GA 
results in smaller pitch transients during the maneuver. 
The crosswind gust capability of the GA effector suite 
was about the same as the original; the GA suite could 
tolerate 29.5 ft/s crossvvind gust, while the original 
suite could withstand a 28 ft/s gust. From this analysis, 
the phase I GA appears to have found a good solution. 
Further testing is required to evaluate the phase II GA 
solution. 


Concluding Remarks 

This paper summarizes several years of research 
supported by the Aircraft Morphing Program at NASA 
Langley Research Center. The paper emphasizes the 
use of MDO techniques applied to the control system 
design process for novel aircraft configurations. One 
such novel aircraft concept is the Lockheed- Martin ICE 
configuration with distributed shape-change devices to 
generate control moments. Starting with a CFD model 
of the ICH configuration, this paper shows how the 
control moments were estimated, promising effector 
locations were proposed, optimal subsets of these 
proposed locations were chosen and the control system 
was simulated and tested. 


Simulation of the Control System 

Several effector suites defined with the MATLAB- 
based tool have been applied to the ICE vehicle in a 
simulation and used in a stability augmentation and 
control system design (Ref. 9). For the present study, 
the control system deployed the effectors in a 


Obviously, the advanced flow control research is not 
complete. The concepts must be tested in the wind 
tunnel, and more effective shape-change dev ices must 
be sought. Moreover, both CFD predictions and flight 
control simulations need to be improved and validated. 
These improved prediction capabilities may influence 
the fitness function used for optimization as well as the 
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implementation details used for automatic 
differentiation. 

While acknowledging that flow control research is in 
its infancy, this paper demonstrates useful MDO 
techniques that will be available to the control 
designers of the future. The automatic differentiation 
techniques demonstrated with the ICE model and 
PMARC code are easily adaptable to other CFD 
models and codes. Likewise, the genetic algorithms 
developed herein can be used with improved control- 
effectiveness measures to find the best locations for a 
wide variety of shape-change effectors. Thus, 
simulation techniques such as CFD, as well as MDO 
techniques such as GA and automatic differentiation, 
empower engineers to explore revolutionary control 
concepts for aircraft of the future. 
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Fig. 1 Control system design process. 
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Fig. 4 One of the 34 candidate effector arrays on the right wing upper surface. 
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Fig. 5 Control authority plots for candidate effector arrays: 7 regions. 
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(a) Population extremes 




(c) Global best 

Fig. 6 Genetic algorithm convergence history. 



(a) Upper surface arrays 



(b) Lower surface arrays 
Fig. 7 Best effector suite found by phase I GA. 
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Fig. 8 Genetic algorithm conv ergence history. 
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(a) Upper surface devices 



Fig. 9 Best set of dev ices found by phase II GA. 
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(b) Pitch moment time histories 



(c) Yaw moment time histories 


Fig. 10 Simulation time histories generated in response 
to a +/- 20-deg bank angle doublet command. 

( ideal — Phase I G A Ref. 9) 


(a) Roll moment time histories 
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